/*
 * Copyright (C) 2014 Nathan Templon
 *
 * This program is free software: you can redistribute it and/or modify
 * it under the terms of the GNU General Public License as published by
 * the Free Software Foundation, either version 3 of the License, or
 * (at your option) any later version.
 *
 * This program is distributed in the hope that it will be useful,
 * but WITHOUT ANY WARRANTY; without even the implied warranty of
 * MERCHANTABILITY or FITNESS FOR A PARTICULAR PURPOSE.  See the
 * GNU General Public License for more details.
 *
 * You should have received a copy of the GNU General Public License
 * along with this program.  If not, see <http://www.gnu.org/licenses/>.
 */

package mathlib;

import function.Function.Domain;
import function.Function.FuncPoint;
import function.Function.Range;
import function.Function.Variable;
import function.SingleVariableTableFunction;
import geometry.Circle;
import geometry.PlaneArea;
import geometry.angle.Angle;
import geometry.angle.Angle.AngleType;
import geometry.angle.Angle.TrigFunction;
import graph.Plot2D;

import java.awt.Color;
import java.util.ArrayList;
import java.util.Collection;

import javax.swing.JFrame;
import javax.swing.JOptionPane;

import log.Logger;
import propulsion.propellant.LiquidPropellant;
import propulsion.rocket.LiquidRocketEngine;
import propulsion.rocket.LiquidRocketEngine.LiquidEngineState;
import propulsion.rocket.LiquidRocketEngine.MassFlowRate;
import propulsion.rocket.component.CombustionChamber;
import propulsion.rocket.component.CylindricalCombustionChamber;
import propulsion.rocket.component.Nozzle;
import propulsion.rocket.component.Throat;
import propulsion.thermo.Gas;
import propulsion.thermo.Gas.HeatRatio;
import propulsion.thermo.Gas.MolarMass;
import propulsion.thermo.IdealGas;
import space.CelestialBody;
import space.LambertSolver;
import space.LambertSolver.LambertSolution;
import space.Orbit;
import space.Orbit.OrbitBuilder;
import space.Orbit.OrbitBuilder.OrbitDefinitionType;
import space.OrbitPosition;
import vector.Vector3D;

/**
 * @author Nathan Templon
 */
public class MathLib {

    public static void main(String[] args) {
        Logger logger = Logger.getDefaultLogger();
        JFrame frame = null;
        
        try {
        
        IdealGas keroseneExhaust = new IdealGas( new Gas.HeatRatio(1.219), new Gas.MolarMass(21.65) );
        IdealGas air = new IdealGas( new HeatRatio(1.4), new MolarMass(28.9));
        Circle throatSection = new Circle( new Circle.Diameter( 0.0254 * (1.5 / 4.0) ) );
        PlaneArea.Area throatArea = throatSection.area();
        Throat throat = new Throat(throatArea);
        Nozzle nozzle = new Nozzle(throatArea, new Nozzle.AreaRatio( 2 ));
        CombustionChamber chamber = new CylindricalCombustionChamber(new Gas.GasState(new Gas.Temperature(3540), new Gas.Pressure(1)));
        
        LiquidRocketEngine redPanda = new LiquidRocketEngine(nozzle, chamber, throat);
        redPanda.setAmbientPres( new Gas.Pressure(101325) );
        
        LiquidPropellant propellant = new LiquidPropellant(keroseneExhaust, new Gas.Temperature(3540));
        LiquidEngineState state = new LiquidEngineState(new MassFlowRate(.05), propellant);
        
//        redPanda.setState(state);
//        
//        String output = "";
//        output += redPanda.thrust() / 4.448 + " Pounds";
//        output += "\n\n" + redPanda.chamberPres().value() / 101325 + " atm chamber";
//        output += "\n\n" + redPanda.exitPres().value() / 101325 + " atm exit";
//        output += "\n\n" + redPanda.isp() + " seconds";
//        
//        popup(output);
        
//        ObjectDB dat = ObjectDB.getGlobalDB();
//        try {
//            dat.addDataFromDir( Settings.defaultDirectory() + "\\Data");
//            dat.buildObjects();
//        }
//        catch (Exception e) {
//            logger.log(e);
//        }
//        
//        
//        Unit fps = (Unit)dat.getObject("UNIT_FEET_PER_SECOND");
//        Unit mps = (Unit)dat.getObject("UNIT_METER_PER_SECOND");
//        Unit mps2 = (Unit)dat.getObject("UNIT_METER_PER_SECOND_PER_SECOND");
//        Unit meter = (Unit)dat.getObject("UNIT_METER");
//        Unit feet = (Unit)dat.getObject("UNIT_FOOT");
//        Unit ft2s = (Unit)dat.getObject("UNIT_FEET_SQUARED_PER_SECOND");
//        Unit m2s = (Unit)dat.getObject("UNIT_METER_SQUARED_PER_SECOND");
//        Unit hertz = (Unit)dat.getObject("UNIT_HERTZ");
//        
//        UnitQuantity x = new UnitQuantity(5, fps);
//        UnitQuantity y = new UnitQuantity(3, feet);
//        UnitQuantity v = new UnitQuantity(2.5, feet);
//        UnitQuantity w = new UnitQuantity(1.75, meter);
        
//        Unit newUnit = feet.times(fps);
//        TreeMap<UnitType, Integer> powMap = hertz.powerMap();
//        String output2 = "";
//        int typeIndex = 0;
//        for(UnitType type : powMap.keySet()) {
//            output2 += type.name() + "^" + powMap.get(type);
//            if (typeIndex < powMap.keySet().size() - 1) {
//                output2 += "\n";
//            }
//            typeIndex++;
//        }
//        
//        TreeMap<UnitType, Integer> powMap2 = (x.dividedBy(y).unit()).powerMap();
//        String output3 = "";
//        int typeIndex2 = 0;
//        for(UnitType type : powMap2.keySet()) {
//            output3 += type.name() + "^" + powMap2.get(type);
//            if (typeIndex2 < powMap2.keySet().size() - 1) {
//                output3 += "\n";
//            }
//            typeIndex2++;
//        }
        
        Variable var = new Variable("x");
        SingleVariableTableFunction func = new SingleVariableTableFunction(var);
        Collection<FuncPoint> funcPoints = new ArrayList<>();
        for(int index = -100; index <= 500; index++) {
            funcPoints.add(new FuncPoint(index / 10.0, index * index / 100.0));
        }
        func.addPoints(funcPoints);
        var.setValue(4.9);
        
        Collection<FuncPoint> points = new ArrayList<>();
        

        
        for(double index = func.domain().min(); index < func.domain().max(); index += 0.1) {
            var.setValue(index);
            double funcValue = func.getValue();
            points.add(new FuncPoint(index, funcValue));
        }
        
        
        CelestialBody earth = new CelestialBody("Earth", 5.9726E24);
        CelestialBody sun = new CelestialBody("Sun", 1988500E24);
        CelestialBody jupiter = new CelestialBody("Jupiter", 1.8981E27);
        
        OrbitBuilder builder = new OrbitBuilder(OrbitDefinitionType.SCALAR_ELEMENT);
        builder.setCentralBody(earth);
        builder.setSMA(6378.0 * 1.1);
        builder.setEcc(0);
        builder.setInc(new Angle(0));
        builder.setRightAscOfAscNode(new Angle(0));
        builder.setArgOfPeri(new Angle(0));
        Orbit parkingOrbit = builder.build();
        OrbitPosition pos = new OrbitPosition(parkingOrbit, new Angle(0));
        
        
        Collection<FuncPoint> bottomPoints = new ArrayList<>();
        for (double f = 180; f < 360; f+= 0.01) {
        	Angle trueAnom = new Angle(f, AngleType.DEGREES);
        	pos.setTrueAnomaly(trueAnom);
        	Vector3D position = pos.radiusVector();
        	double xPos = position.getComponent(1);
        	double yPos = position.getComponent(2);
        	bottomPoints.add(new FuncPoint(xPos, yPos));
        }
        Variable bottomVar = new Variable("X");
        SingleVariableTableFunction bottomFunc = new SingleVariableTableFunction(bottomVar);
        bottomFunc.addPoints(bottomPoints);
        Collection<FuncPoint> topPoints = new ArrayList<>();
        for (double f = 0; f < 180; f+= 0.01) {
        	Angle trueAnom = new Angle(f, AngleType.DEGREES);
        	pos.setTrueAnomaly(trueAnom);
        	Vector3D position = pos.radiusVector();
        	double xPos = position.getComponent(1);
        	double yPos = position.getComponent(2);
        	topPoints.add(new FuncPoint(xPos, yPos));
        }
        Variable topVar = new Variable("X");
        SingleVariableTableFunction topFunction = new SingleVariableTableFunction(topVar);
        topFunction.addPoints(topPoints);
        
        
        OrbitBuilder surfaceBuilder = new OrbitBuilder(OrbitDefinitionType.SCALAR_ELEMENT);
        surfaceBuilder.setCentralBody(earth);
        surfaceBuilder.setSMA(6378.0);
        surfaceBuilder.setEcc(0);
        surfaceBuilder.setInc(new Angle(0));
        surfaceBuilder.setRightAscOfAscNode(new Angle(0));
        surfaceBuilder.setArgOfPeri(new Angle(0));
        Orbit surfaceOrbit = surfaceBuilder.build();
        OrbitPosition surfacePos = new OrbitPosition(surfaceOrbit, new Angle(0));
        
        Collection<FuncPoint> bottomSurfPoints = new ArrayList<>();
        for (double f = 180; f < 360; f+= 0.01) {
        	Angle trueAnom = new Angle(f, AngleType.DEGREES);
        	surfacePos.setTrueAnomaly(trueAnom);
        	Vector3D position = surfacePos.radiusVector();
        	double xPos = position.getComponent(1);
        	double yPos = position.getComponent(2);
        	bottomSurfPoints.add(new FuncPoint(xPos, yPos));
        }
        Variable bottomSurfVar = new Variable("X");
        SingleVariableTableFunction bottomSurfFunc = new SingleVariableTableFunction(bottomSurfVar);
        bottomSurfFunc.addPoints(bottomSurfPoints);
        Collection<FuncPoint> topSurfPoints = new ArrayList<>();
        for (double f = 0; f < 180; f+= 0.01) {
        	Angle trueAnom = new Angle(f, AngleType.DEGREES);
        	surfacePos.setTrueAnomaly(trueAnom);
        	Vector3D position = surfacePos.radiusVector();
        	double xPos = position.getComponent(1);
        	double yPos = position.getComponent(2);
        	topSurfPoints.add(new FuncPoint(xPos, yPos));
        }
        Variable topSurfVar = new Variable("X");
        SingleVariableTableFunction topSurfFunction = new SingleVariableTableFunction(topSurfVar);
        topSurfFunction.addPoints(topSurfPoints);
        
        
        Vector3D hyperbolicRadius = new Vector3D(0, -1.1 * 6378, 0);
        Vector3D hyperbolicVelocity = new Vector3D(13, 0, 0);
        OrbitPosition hyperbola = new OrbitPosition(earth, hyperbolicRadius, hyperbolicVelocity);
        Collection<FuncPoint> hyperPoints = new ArrayList<>();
        for (double f = -110; f < 110; f+= 0.01) {
        	Angle trueAnom = new Angle(f, AngleType.DEGREES);
        	hyperbola.setTrueAnomaly(trueAnom);
        	Vector3D position = hyperbola.radiusVector();
        	double xPos = position.getComponent(1);
        	double yPos = position.getComponent(2);
        	hyperPoints.add(new FuncPoint(xPos, yPos));
        }
        Variable hyperVar = new Variable("X");
        SingleVariableTableFunction hyperFunc = new SingleVariableTableFunction(hyperVar);
        hyperFunc.addPoints(hyperPoints);
        
//        Collection<FuncPoint> pandaPoints = new ArrayList<>();
//        for(double ratio = -100; ratio <= 100; ratio += 0.1) {
//            pandaPoints.add(new FuncPoint(ratio, ratio * ratio * 0.5));
//        }
//        Variable pandaVar = new Variable("SemiMajorAxis");
//        SingleVariableTableFunction pandaFunc = new SingleVariableTableFunction(pandaVar);
//        pandaFunc.addPoints(pandaPoints);
        
        
        Plot2D graph = new Plot2D();
        graph.setNumFuncPoints(500);
        graph.setPreferredSize(graph.getPreferredSize());
        
//        double pandaDomainSpan = pandaFunc.domain().max() - pandaFunc.domain().min();
//        double pandaRangeSpan = pandaFunc.range().max() - pandaFunc.range().min();
//        graph.setDomain(new Domain(pandaFunc.domain().min() - pandaDomainSpan * 0.1, pandaFunc.domain().max() + pandaDomainSpan * 0.1));
//        graph.setRange(new Range(pandaFunc.range().min() - pandaRangeSpan * 0.1, pandaFunc.range().max() + pandaRangeSpan * 0.1));
        graph.setDomain(new Domain(-50000, 50000));
        graph.setRange(new Range(-50000, 50000));
        graph.setXTicksPerGridline(5);
        graph.setYTicksPerGridline(5);
        graph.setPointRadius(5);
        graph.setPointColor(Color.BLACK);
        graph.setLineColor(Color.BLACK);
        
        graph.setXAxisTitle("X");
        graph.setYAxisTitle("Y");
        
//        graph.addFunction(pandaFunc, Plot2D.FUNC_TYPE.LINE);
//        graph.addFunction(pandaFunc, Plot2D.FUNC_TYPE.POINT);
        graph.addFunction(bottomFunc, Plot2D.FUNC_TYPE.LINE);
        graph.addFunction(topFunction, Plot2D.FUNC_TYPE.LINE);
        graph.addFunction(bottomSurfFunc, Plot2D.FUNC_TYPE.LINE);
        graph.addFunction(topSurfFunction, Plot2D.FUNC_TYPE.LINE);
        graph.addFunction(hyperFunc, Plot2D.FUNC_TYPE.LINE);
        graph.updateFunctions();
        
        frame = new JFrame("Test of Plot2D");
        frame.setDefaultCloseOperation(JFrame.EXIT_ON_CLOSE);
        frame.getContentPane().add(graph);
        frame.pack();
//        frame.setVisible(true);
        
        graph.square();
        
        
        final double AU = 149597871.0;
        final double YEAR_S = 1.0 / (60 * 60 * 24 * 365);
        
        OrbitBuilder build = new OrbitBuilder(OrbitDefinitionType.SCALAR_ELEMENT);
        build.setCentralBody(sun);
        build.setSMA(1.0 * AU);
        build.setEcc(0);
        build.setInc(new Angle(0));
        build.setRightAscOfAscNode(new Angle(0));
        build.setArgOfPeri(new Angle(0));
        Orbit earthsOrbit = build.build();
        OrbitPosition helioStartPos = new OrbitPosition(earthsOrbit, new Angle(0));
        
        build.setCentralBody(sun);
        build.setSMA(5.2 * AU);
        build.setEcc(0);
        build.setInc(new Angle(0));
        build.setRightAscOfAscNode(new Angle(0));
        build.setArgOfPeri(new Angle(0));
        Orbit jupitersOrbit = build.build();
        OrbitPosition helioEndPos = new OrbitPosition(jupitersOrbit, new Angle(Math.PI));
        
        LambertSolver transferSolver = new LambertSolver(helioStartPos, helioEndPos);
        LambertSolution transferSolution = transferSolver.solution(transferSolver.aMin(), true);
        
        final Vector3D k = new Vector3D(0, 0, 1);
        
        Vector3D position = transferSolution.finalPosition().radiusVector();
        Vector3D va = transferSolution.secondVelocity();
        Vector3D jupVel = helioEndPos.velocityVector();
        Vector3D jupRadius = helioEndPos.radiusVector();
        
        Vector3D vInfVec = (Vector3D)va.add(jupVel.scalarMultiply(-1.0));
        double vInf = vInfVec.magnitude();
        
        final double rJup = 71370.0;
//        for (int i = 208000; i < 209000; i++) {
        	double rp = 20.8934 * rJup;
        	double sinDeltaBy2 = sinDeltaBy2(vInf, rp, jupiter);
        	Angle delta = (new Angle(sinDeltaBy2, TrigFunction.SINE)).scalarMultiply(2.0);
        	// Directing the s/c to fly under Jupiter and be deflected upwards
        	delta = new Angle(2.0 * Math.PI - delta.measure());
        	Vector3D hHat = (Vector3D)position.cross(va).direction();
        	Vector3D vHat = (Vector3D)va.direction();
        
        	Vector3D vInfOut = (Vector3D) vHat.scalarMultiply(vInf * delta.cos()).add(hHat.scalarMultiply(vInf * delta.sin()));
        	Vector3D vFinal = (Vector3D) vInfOut.scalarMultiply(-1.0).add(jupVel);
        
        	OrbitPosition newOrbit = new OrbitPosition(sun, jupRadius, vFinal);
        	double perRatio = newOrbit.getOrbit().period() / jupitersOrbit.period();
        	if (1.4999 <= perRatio && perRatio <= 1.5001) {
        		System.out.println(perRatio);
        	}
//        }
        	
       System.out.println("V_FINAL: " + vFinal);
       System.out.println("Inclination = " + newOrbit.getInc().measure(AngleType.DEGREES) + " Degrees");
       System.out.println("Semi Major Axis = " + (newOrbit.getSMA() / AU) + " AU");
       System.out.println("Aphelion = " + (newOrbit.getOrbit().getApo() / AU) + " AU");
       System.out.println("RightAscOfAscNode = " + newOrbit.getRightAscOfAscNode().measure(AngleType.DEGREES) + " Degrees");
       System.out.println("Argument of Peri = " + newOrbit.getArgPeri().measure(AngleType.DEGREES) + " Degrees");
       System.out.println("DELTA = " + delta.measure(AngleType.DEGREES));
        
        Vector3D helioVelVel_2 = newOrbit.radiusVector();
        
        }
        catch(Exception e) {
            if (frame != null) {
                frame.dispose();
            }
            logger.log(e);
            JOptionPane.showMessageDialog(null, e);
        }
        
        logger.close();
    }
    
    public static void popup(Object o) {
        JOptionPane.showMessageDialog(null, o);
    }

    
    public static double sinDeltaBy2(double vInf, double rp, CelestialBody body) {
    	return 1 / ( 1 + ( (vInf * vInf * rp) / body.gravParam() ) );
    }
    
    
}
